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Behavior analogous to a second order phase transition is observed for the homogeneous cooling 
state of a heavy impurity particle in a granular fluid. The order parameter is the ratio of impurity 
mean square velocity to that of the fluid, with a conjugate field h proportional to the mass ratio. 
A parameter /3, measuring the fluid cooling rate relative to the impurity-fluid collision rate, is the 
analogue of the inverse temperature. For 13 <1 the fluid is "normal" with = at /i = 0, as in the 
case of a system with elastic collisions. For /3 > f an "ordered" state with 7^ occurs at /i = 0, 
representing an extreme breakdown of equipartition. Critical slowing and qualitative changes in the 
velocity distribution function for the impurity particle near the transition are noted. 

PACS number(s): 45.70.-n, 45.70.Mg, 05.20.Dd 



The simplest statistical mechanical model for an acti- 
vated granular fluid is a system of smooth, inelastic hard 
spheres. The inelasticity is specified in terms of a resti- 
tution coefficient a < 1. Qualitative differences between 
the cases a — 1 (elastic spheres) and a ^ 1 have been 
demonstrated in many ways using both theoretical and 
simulation methods . Among these is the replacement 
of the Gibbs state for an isolated uniform system by a ho- 
mogeneous cooling state (HCS). The HCS is inherently 
time dependent due to the loss of energy in each binary 
collision. It is postulated that the time dependence of 
this ensemble occurs entirely through a scaling of the ve- 
locity by its root mean square ( "thermal" ) velocity vt (t) 
and the associated normalization In this case, the 
dynamics (Liouville equation or kinetic equation) can be 
transformed to a stationary state form using suitable di- 
mensionless variables. The time dependence of VT{t) also 
follows from this scaling, so that only a time indepen- 
dent dimensionless cooling rate must be determined self- 
consistently. Similar results apply for a multi-component 
system, although the analysis is complicated by the ex- 
istence of many scaling velocities The HCS results 
when the cooling rates for all mean square velocities are 
the same. This condition then determines all scaling ve- 
locities in terms of that for one species. It is the analogue 
of energy equipartition for a system with a — I, but de- 
viations from equipartition occur for all a 7^ 1. These 
deviations can be weak or strong depending on the me- 
chanical differences between the species and the degree 
of dissipation. 

Many of the new effects to be expected in multi- 
component systems can be illustrated by the simplest 
case of a single impurity particle in a one component 
system [Q. The fiuid (impurity) particle diameter, mass, 
and restitution coefficient are denoted by a (ctq) , m (mo), 
and a (ag). There are two mean square velocities and 
associated cooling rates for the fluid and impurity parti- 
cles. 
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The brackets denote an average over the HCS, and the 
subscript denotes properties of the impurity particle. 
The cooling rates for a moderately dense fluid can be 
calculated to good accuracy from the Boltzmann-Enskog 
and Lorentz-Enskog kinetic equations for the fluid and 
impurity reduced distribution functions These rates 
then depend on the velocity distributions for the fluid 
and impurity particle in the HCS and can be estimated to 
good approximation using Maxwellians (see justification 
below) with the results 
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where n is the fluid density, g is the pair correlation 
function for two fluid particles, and na^gvT is a char- 
acteristic collision rate for the fluid particles. Similarly, 
ly* — v/ncj'^gvT — (16-\/7r/3)ft. (50/5) (o'/cr)^ is a corre- 
sponding dimensionless average impurity-fluid particle 
collision rate, where a — {a + aa) /2 and go is the cor- 
responding fluid-impurity particle pair correlation func- 
tion. The dependence on the mass ratio occurs entirely 
through h = to(1+q;o) /2(m-|-mo). The time independent 
ratio of mean square velocities. 
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is determined by the condition that the cooling rates be 
equal, 
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The constant /3 is a measure of the fluid cooling rate 
relative to the fluid-impurity collision rate. 
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The elastic collisions limit is given hy f3 — 0, for which 
the solution is (p = h/ {1 — h) — m/mQ. This is the 
required result from the equipartition theorem. If only 
the fluid particle collisions are elastic {(3 = 0, ao 7^ 1); 
a recent exact result of Martin and Piasecki, (j> = m(l + 
ao)/ [2mo + m (1 — ao)] is recovered ||]. More generally, 
Eq. ^ determines h) for the entire physical domain 
of h and /3. The solution is found to have qualitatively 
different behavior at small h for different values of f3. The 
asymptotic behavior for — s- is found to be 

({l-P)-'h, /3<1 

<^ V2h, (3=1 (7) 

[ ~l + 2(3^{(3^ -l)-'^h, 13 > I. 

Thus (p vanishes at /i = for /3 < 1 but is finite at /3 > 1. 
The condition /3 < 1 is "normal" in the sense that 
scales as the mass ratio just as for systems with elas- 
tic collisions and the consequent equipartition of kinetic 
energy. This normal state has been studied recently for 
/i <C 1, where the kinetic theory implies a correspond- 
ing Fokker-Planck description Q]. The condition (3 > 1 
corresponds to an extreme violation of equipartition such 
that the mean square velocities of the fluid and impurity 
particles remain comparable even though their mass ratio 
is zero. Figure ^ shows the more complete dependence of 
(j>{(3, h) on /i for /? = 0.9, 1, and 1.1, confirming the above 
asymptotic analysis. 
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FIG. 2. Inverse susceptibility expansion coefficient 

(ah), and heat capacity (Ch) as functions of /3 for h = 10~^ 
(dotted line), h = 10~^ (dashed line), and /i = (solid line). 



FIG. 1. Ratio of mean square velocities, </!>, as a function of 
the mass ratio parameter h for /3 = 0.9, 1, and 1.1. 



This behavior at small h is analogous to a thermody- 
namic transition at /3 = 1 between two different states 
of the system characterized by the order parameter 0, 
conjugate field h, and inverse "temperature" (3. To de- 
velop this analogy a Helmholtz free energy F {f3, (j>) is 
obtained from the "equation of state" (||) by integrating 
the definition h = OF {/3, (f) /90 [with the boundary con- 
dition F{(3,0) = 0, for simplicity]. Next, the Gibbs free 
energy is constructed from F by the Legendre transfor- 
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mation <I> (/?, h) = F — h(j>. The first and second deriva- 
tives of $ {(3, h) provide the order parameter (j), "entropy" 
S, "susceptibihty" x^, "expansion coefheient" ah, and 
"heat capacity" Ch- The results are 
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The behavior of these functions near h = foUows di- 
rectly from (|^). For example, the response functions X/3j 
ah, and at ft- = 0, 10^^^ a^ji 10~^ are shown in Fig. 
Ij. In the limit h 0, the entropy and order param- 
eter are continuous a,t P — 1, so the divergence of xp 
and the discontinuity of ah at /3 = 1 identifies this as a 
second order phase transition. Near the critical region 
{h <C 1, 1/3 — 1| <C 1), the free energy adopts the Landau- 
like form $(/3, /i) « i(l — /3)02 + i03 — h(j), which yields 
h ^ {1 — f3)(f) + ^(j)^ . As a consequence, the free energy 
and the equation of state in the critical region satisfy 
the scaling relations $(A(/3 - l),X"^h) = A^$(/3 - l,h) 
and 0(A(/3 - 1), A'^ft) = A''-'' (/>(/? - 1, h), with a = 2 and 

6 = 3. These scaling relations suffice to determine the 
critical exponents |^ ^ = a/{b—a) = 2, /3 = b—a = 1, and 

7 — 2a — b — 1, while the critical exponent a — 2 — 6= —1 
is negative, indicating that Ch is continuous at the criti- 
cal point. 

If the ratio between the initial mean square velocities 
of the fluid and impurity particles is not that given by the 
solution to (||), there is an evolution to the HCS described 
by (^ which can be written in the Ginzburg-Landau form 



9$ (/3,/i; 



(13) 



where a new time scale ds — na'^gvT(t)dt has been intro- 
duced. Here, $ (/3, h; (j)) is a variational free energy given 
by Eq. (|^) with the order parameter </> considered as an 
independent variable, and the kinetic coefficient /i((/>) is 
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The stationary solution occurs for 9$ (/?, h; (f) /90 = 0, 
which is just Eq. (||). For states near the HCS the equa- 
tion can be linearized and a characteristic response time 
r identified according to 
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where i^hcs denotes the solution to (|^). In the elastic 
limit r is just the equilibration time for the impurity 
particle to attain a mean kinetic energy equal to that of 
the fluid particles. Similarly, for inelastic collisions it is 
the time (in units of the fluid mean free time) for the 
impurity particle to reach a cooling rate equal to that of 
the fluid. This characteristic time is a smooth function 
of h and /? except in the limit /i — > where r diverges 
ai P ~ 1. This critical slowing follows directly from the 
fact that r (X X/3- Otherwise, the relaxation times away 
from P = 1 are finite and comparable for the normal and 
ordered states. Figure ^ shows the dependence of the 
inverse time t*~^ = {v*t)~^ on ft for /3 = 0.9, 1, and 1.1. 
The quantity r* measures the relaxation time for the 
impurity particle in terms of the number of impurity- 
fluid particle collisions. 
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FIG. 3. Inverse characteristic time r* = (i^'t) ^ as a 
function of the mass ratio parameter h for l3 = 0.9, 1, and 
1.1. 

The above analysis is based on an approximate evalu- 
ation of the cooling rates using Maxwellians for the im- 
purity and fluid particle velocity distributions. This is 
known to be accurate for the fluid particles but the 
new features in the critical domain and in the ordered 
phase /? > 1 require closer inspection of the impurity 
particle distribution. Based on the Lorentz-Enskog ki- 
netic equation, it is possible to determine the impurity 
particle velocity distribution for the HCS asymptotically 
close to /? ^ 1 and ft ^ with the result 
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where v* = Vq/vto and C is a normalization constant. 
The asymptotic equation of state is determined from 
/ dv* {v*^ - I) F(v*) = 0. The equation of state is 
similar to that of Eq. with only minor quantitative 
changes in some coefficients. This justifies the preced- 
ing analysis based on the Maxwellian approximation to 
fo- Nevertheless, the distribution function itself shows 
new qualitative differences. As expected, for (3 < 1, (p 
h/ {1 — (3) and the Maxwellian is recovered. For f3 > 1, 
(f) 10 (/3 — 1) /3 and the velocity distribution becomes 

sharp at v* = VW2, ^(v*) (67r)-M (w* - ^3/2) • 

Furthermore, at the critical point (/? = 1), ^ — > VlOA/i 
and F(v*)) ^ Cexp(-AD*'*), where A ~ 0.24. This 
crossover from Gaussian to exponential quartic to delta 
distributions is illustrated in Fig. ^ for /i = 10~^ and 
(i = 0.9, 1.0, and 1.1. 
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FIG. 4. Velocity distribution function of the impurity par- 
ticle, F, relative to the Maxwellian, Fm-, for h = 10~^ and 
13 = 0.9 (dotted line), (3 = 1 (dashed line), and /3 = 1.1 (solid 
line) . 

The origin of this transition can be traced to the fact 
that the cooling of the fluid is due to collisions that de- 
pend on the restitution coefficient a but not the mass 
ratio (or h), while cooling of the impurity depends on 
the mass ratio but not a. The equality of cooling rates 
is enforced by adjusting the ratio of the mean square 
velocities, (f>. In general, for fixed mechanical proper- 
ties, this is possible for small mass ratio only if (j) is non 
zero (the ordered state); the cooling rate for the impu- 
rity has an explicit proportionality to h which must be 
countered by (p. This implies that the mean kinetic en- 
ergy of the impurity particle diverges relative to that of 
the fluid particles. However, if this explicit decrease of 
the impurity cooling rate with h is tempered by a large 
size ratio and/or weak fluid particle dissipation then the 
normal behavior with </) = is restored. The parame- 
ter [3 quantifies the meaning of "tempered" . To give an 



explicit example, consider the case a = ag = 0.95 and 
to/ mo = 10^^. Also, for simplicity consider low den- 
sity so that g = gQ 1. Then j3 — 1.79 (cr/a") and 
the system is close to the ordered state for cr = cto but 
close to the normal state for ctq = 2(T. It is possible 
to study the crossover through the critical domain using 
Direct Monte Carlo simulation of the Boltzmann-Enskog 
and Lorentz-Enskog kinetic equations, and by molecular 
dynamics simulation. In fact both Monte Carlo molec- 
ular dynamics simulations for the domain (3 < 1 and 
m/mQ w 10~^ have been performed to confirm Brownian 
motion in the normal state . Extension of these studies 
to the critical and ordered states should be straightfor- 
ward. 

For experimental purposes it is perhaps more relevant 
to study these critical properties as refiected in the diffu- 
sion coefficient D (or mean square displacement) for the 
impurity particle. The dependence of D on /3 and h can 
be calculated from the Lorentz-Enskog kinetic equation 
in the same approximation as the cooling rates above. 
The studies of the Fokker-Planck limit (/i ^ 0) in Q sug- 
gest a divergence of Z) as /? ^ 1 although that analysis 
is strictly limited to /3 < 1. A more complete description 
of D for all /3, h will be given elsewhere |^ . 
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